Correlated electron-nuclear dynamics: Exact factorization of the molecular wavefunction 
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It was recently shown fj] that the complete wavefunction for a system of electrons and nuclei 
evolving in a time-dependent external potential can be exactly factorized into an electronic wavefunc- 
tion and a nuclear wavefunction. The concepts of an exact time-dependent potential energy surface 
(TDPES) and exact time-dependent vector potential emerge naturally from the formalism. Here we 
present a detailed description of the formalism, including a full derivation of the equations that the 
electronic and nuclear wavefunctions satisfy. We demonstrate the relationship of this exact factoriza- 
tion to the traditional Born-Oppenheimer expansion. A one-dimensional model of the molecule 
in a laser field shows the usefulness of the exact TDPES in interpreting coupled electron-nuclear dy- 
namics: we show how features of its structure indicate the mechanism of dissociation. We compare 
the exact TDPES with potential energy surfaces from the time-dependent Hartree-approach, and also 
compare traditional Ehrenfest dynamics with Ehrenfest dynamics on the exact TDPES. 



I. INTRODUCTION 



The interplay of nuclear and electronic dynamics in 
the presence of time-dependent external fields leads to 
fascinating phenomena, especially beyond the perturba- 
tive regime, e.g. photo-induced molecular dissociation, 
charge-resonance enhanced ionization, control of elec- 
tron localization, electron-hole migration after photo- 
excitation, to name a few 0-@]. The exact solution of 
the time-dependent Schrodinger equation (TDSE) is cur- 
rently out of computational reach except for the very 
simplest of molecules [7], such as H^~, so usually ap- 
proximate methods are used. Typically, (but not always, 
see Refs. l8l-[TTll), these methods treat the nuclei clas- 
sically as point charges with electron-nuclear coup ling 
given by Ehrenfest dynamics, or surface-hopping 1131 ; 
a topical a ppl ication is to model photochemical pro- 
cesses EUTlll] , for example, in solar cells, to study 
the (field-free) dynamics ensuing after an initial elec- 
tronic excitation. Indeed several examples have shown 
that the predicted electron-hole migration can depend 
critically on the description of the nuclear motion and 
how it is correlated with the electronic dynamics (see 
Ref. (a, la] and references within). Apart from enabling 
calculations on more than the simplest systems possi- 
ble, these methods provide much intuition, in particu- 
lar through the central concept of the potential energy 
surface (PES). Indeed, the very idea itself of surface- 
hopping would not exist without the notion of a land- 
scape of coupled PESs. Dressed molecular potentials 
such as light-induced molecular potentials (LIMPS) fl5ll 
have proved valuable in understanding processes such 
as bond-softening, stabilization against dissociation, etc. 
where the laser field induces avoided crossings between 
PESs. Approximate time-dependent potential energy 
surfaces (TDPES) were introduced by Kono lfl6ll as in- 
stantaneous eigenvalues of the electronic Hamiltonian, 



and have proven extremely useful in the interpretation 
of system-field phenomena, as have the quasi-static or 
phase-adiabatic PES's used recently to interpret electron 
localization in dissociative ionization [17]. Recent work 
of Cederbaum [18] introduced a TDPES in a different 
way, by generalizing the Born-Oppenheimer approxi- 
mation to include time-dependent external potentials. 
In short, the PES is perhaps the most central concept in 
our understanding of molecular motion. 

In a recent Letter we showed that an exact TDPES 
may be defined, via a rigorous separation of electronic 
and nuclear motion by introducing an exact factorization 
of the full electron-nuclear wavefunction. The idea of an 
exact factorization was first introduced by Hunter fl9ll 
for the static case. He also deduced the exact equa- 
tion of motion for the nuclear factor. The equation of 
motion for the electronic wavefunction was first given 
by Gidopoulos and Gross j20ll for the time-independent 
case. Both in the static and in the time-dependent case 
the factorization leads to an exact definition of the PES, 
and also of the Berry vector potential. What is par- 
ticularly interesting about the vector potential is that 
Berry-Pancharatnam phases felll are usually interpreted 
as arising from some approximation where a system is de- 
coupled from "the rest of the world", thereby making 
the system Hamiltonian dependent on some "environ- 
mental" parameters. For example, in the static BO ap- 
proximation, the electronic Hamiltonian depends para- 
metrically on nuclear positions, and when the molecular 
wavefunction is approximated by a single product of a 
nuclear wavefunction and an eigenstate of the electronic 
Hamiltonian, the equation of motion for the former con- 
tains a Berry vector potential. The question whether the 
BO Berry phase survives in the exact treatment was first 
discussed in Ref. ll20Tl for the static case and in Ref. Q]] 
for the time-dependent case. 

In the present paper, we provide the detailed deriva- 
tion of the formalism of Ref. Q]] (Section analyse 



features of the exact electron-nuclear coupling terms in 
general (Section [TTT|> . including their relationship to cou- 
plings in the traditional Born-Oppenheimer expansion, 
and then study the TDPES for the specific case of a 
model molecule in an oscillating electric field (Sec- 
tion llV) . The remainder of this introduction serves to 
set up the problem at hand, and to remind the reader 
of the Born-Oppenheimer treatment of the electron- 
nuclear system. 



A. The Hamiltonian 



The particle-particle Coulomb interactions have the 
form: 
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1,3 = 1 



(7) 
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In this section we establish notation and define the 
Hamiltonian for the combined system of electrons and 
nuclei. The coordinates of the N e electrons are collec- 
tively denoted by r s where r = {rj} and s = {sj},j = 
l...N e , represent electronic spatial and spin coordinates, 
respectively The N n nuclei have masses M 1 ...Mj^ n and 
charges Z\...Zff n and coordinates collectively denoted 
by R (7 where R = {R a } and a_ = a = l...N„, rep- 

resent nuclear spatial and spin coordinates, respectively. 
Furthermore, we consider the system is under the influ- 
ence of some time-dependent external scalar field. The 
system is described, non-relativistically by the Hamilto- 
nian 

H = H BO (r, R) + f c e xt (r, t) + f n (R_) + (E, t), (1) 

where Hbo(l, R) is the familiar Born-Oppenheimer 
electronic Hamiltonian, 
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The quantum mechanical equation of motion of such 
a system is given by the TDSE: 



#*(r s,Ra,() = !3 ( $(rs,Rff,i) 



(10) 



The full electron-nuclear wavefunction, $(rs,R cr, t), 
that satisfies the TDSE flOt , contains the complete in- 
formation on the system. As discussed in the introduc- 
tion, it can be solved numerically only for very small 
systems of one or two electrons and nuclei and , more- 
over, VP does not give access to PESs, which provide an 
intuitive understanding and interpretation of the cou- 
pled electron-nuclear dynamics. 



H BO =T e (|) + W ee (|) + Wen (r, R) + W nn (RJ . (2) 

The subscripts "e" and "n" refer to electrons and nu- 
clei, respectively, and atomic units are used throughout 

(e 2 = h = m e = 1). Here 



^ = -Ek 2 



(3) 



and 



N " 1 

^ = ~ E 



(4) 



denote the kinetic-energy operators of the electrons and 
nuclei, respectively. All external scalar potentials on the 
system (e.g. electric fields) are represented by 



K'xt = E V ext(IW), 



(5) 



and 



(6) 



B. The Born-Oppenheimer Approximation 

The Born-Oppenheimer (BO) approximation is 
among the most basic approximations in the quantum 
theory of molecules and solids. Consider the case 
when there is no external time-dependence in the 
Hamiltonian. The BO approximation relies on the fact 
that electrons typically move much faster than the 
nuclei; on the timescale of nuclear motion, the electrons 
"instantly" adjust to remain on the instantaneous 
eigenstate. This "adiabatic approximation" allows 
us to visualize a molecule or solid as a set of nuclei 
moving on the PES generated by the electrons in a 
specific electronic eigenstate. The electronic Hamilto- 
nian Hbo(Li S) depends parametrically on the nuclear 
positions, via the electron-nuclear Coulomb interaction. 
That is, the stationary electronic Schrodinger equation 
is solved for each fixed nuclear configuration R cr, 



H BO (L, S)4,(o) 



(is) (11) 



yielding (R (j)-dependent eigenvalues V BO (R a) and 
eigenf unctions 4> 3 R . The total molecular wavefunction, 
^ bo (l s, R <L), is then approximated as a product of the 
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relevant electronic state, ^ (r s), and a nuclear wave- 
function x^/°(IL 2l) satisfying the corresponding BO nu- 
clear Schrodinger equation 

(E^x dbH v « + ^(E^)) 2 + 4>(§g)) xf„ (ag) 

= ^°(l£) (12) 



where 



(13) 



Es(0R CT 



J? BO (r,R) + E 



^R <T 



and 



(14) 



where (..|..|..) r denotes an inner product over all spa- 
tial electronic variables only. The index v of the nuclear 
wave function labels the vibrational/rotational eigen- 
state on the jth PES. The second term on the right of 
Eq. [13] is often referred to as the "BO diagonal correc- 
tion" or "adiabatic correction". The potential energy 
surface e J BO (R q_) is enormously important in molec- 
ular physics and quantum chemistry. It is a central 
tool in the analysis and interpretation of molecular ab- 
sorption and emission spectra, experiments involving 
nuclear motion, mechanisms of dissociation, energy- 
transfer, for example. The nuclear dynamics on a single 
PES (sometimes called "BO dynamics") is obtained by 
using the Hamiltonian on the left of Eq. JT2] | in a time- 
dependent Schrodinger equation for a time-dependent 
nuclear wavefunction x(R a, t). This corresponds to 
approximating the total molecular wavefunction by a 
time-dependent nuclear wavepacket multiplied with a 
static electronic BO state: 



*(i |, E g, t) « x BU (E g, t)K a (r |) 



(15) 



The vector potential .F^(R S.), especially the Berry 
phase associated with it, ^J"^(R a) ■ dR, captures 
the essential features of the behavior of a system with 
conical intersections. Inclusion of the Berry phase can 
significantly shift and re-order the energy eigenvalues 
of molecular roto-vibrational spectra, as well as scat- 
tering cross-sections (although sometimes undetected in 
experiments that measure integrated quantities, due to 
cancellations between paths, see e.g. Refs. j22l - l26ll and 
references within). 

It appears from the above discussion that in the tra- 
ditional treatment of molecules and solids the concepts 
of the PES and the Berry phase arise as a consequence 
of the BO approximation. Some of the most fascinat- 
ing phenomena of condensed-matter physics, like su- 
perconductivity, however, appear in the regime where 



the BO approximation is not valid; likewise typical pho- 
todynamical processes in molecules require going be- 
yond the single-electronic-surface picture. This raises 
the question: If one were to solve the Schrodinger equa- 
tion of the full electron-nuclear Hamiltonian exactly (i.e. 
beyond the BO approximation) do the Berry phase and 
the potential energy surface survive, with a possibly 
modified form, and if so, how and where do they show 
up? What is their relation to the traditional potential en- 
ergy surface and Berry phase in the BO approximation? 
Moreover, many interesting phenomena occur when 
molecules or solids are exposed to time-dependent ex- 
ternal field e.g. lasers. Can one give a precise meaning 
to a time-dependent potential energy surface and a time- 
dependent vector potential? 

Before answering the points raised above, focussing 
on the time-dependent case, we briefly discuss the 
Born-Oppenheimer expansion which solves the full 
TDSE Eq. ( TTUT l exactly for the coupled electron-nuclear 
system. 



C. The Born-Oppenheimer Expansion 

The set of electronic eigenf unctions a (r_ s)} calcu- 
lated from Eq. (TTTb form a complete orthonormal set in 
the electronic space for each fixed R q_ 

E / dv Jk 2d i)4 ad 1) = s n > ( 16 ) 

therefore the total time-dependent wavefunction of the 
system ^(r s, R <L t) can be expanded in that basis: 

00 



Here 



xf° 



(Rg,t)=J2 J «fetfS s (Ei)*(kl»B:£,t) (18) 



are the expansion coefficients which are functions of 
the nuclear degrees of freedom and time. Eq. l fl7l l is 
the so-called BO expansion which is an exact repre- 
sentation of the complete molecular wavefunction due 
to the completeness of {0R CT (r s)}. It applies also to 
fully-time-dependent problems where ^ evolves un- 
der external time-dependent potentials V^t- m prac- 
tice, for numerically feasible calculations, approxima- 
tions are introduced to limit the expansion to a small 
subset of {(f) 3 R (r s)}. By inserting the expansion (17) 

into Eq. ((Tol l, multiplying by tfi-^ (r s) from the left, and 
integrating over the electronic degrees of freedom, equa- 
tions for the expansion coefficients X? (E £> t) are d e ~ 
termined. One obtains: 
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[E^Hv«+^ fe °j 2 + 



E[<^i^tWi^'>-E A 



BO 

kj,a 



BO 



BO _ ,M. 
Xj ~ dt 



(19) 



Here 



2M Q 



(20) 



is the time-dependent scalar potential and is the fcth gen- 
eralized BO potential energy, generalized to account for 
the time-dependent external field (c.f. Eq. {13)). The 
terms 



(21) 



are called the "nonadiabatic couplings", defined by fH- 



J 7 ^® = -i<^|V a 4a> 
^° a (R) =<<& ff |V^ R(T > 



(22) 



II. EXACT FACTORIZATION OF THE 
TIME-DEPENDENT ELECTRON-NUCLEAR 
WAVEFUNCTION 

The BO expansion Eq. ((TTt yields the complete molec- 
ular wavefunction exactly. Instead of having an infinite 
sum of terms involving an infinite set of generalized 
PES's and non-adiabatic couplings, the question arises 
whether it is possible to represent the complete, time- 
dependent, electron-nuclear wavefunction exactly as a 
single product of an electronic wavefunction and a nu- 
clear wavefunction. In this section, we show that the 
answer is yes. We derive formally exact equations of 
motion for each subsystem, out of which emerge rig- 
orous definitions of a time-dependent potential energy 
surface (TDPES) and a time-dependent vector potential. 

Visually, the decomposition is similar in form to the 
single-surface BO approximation, yet it is exact. There is 
no assumption on the time scale of the motions of each 
subsystem, i.e. unlike in the BO approximation, we do 
not solve for the "fast" variables first and then feed it 
into the equation for the "slower" variables. Instead, 
the equations of motion for each subsystem are derived 
together, in a variational approach. The exact decom- 
position, contrary to the BO separation, accounts for the 
full correlation between the two subsystems, regardless 
of the mass and energy of the nuclear subsystem. In the 
following we formalize the idea as a theorem which we 
then prove. We discuss in detail the implications of this 
exact decomposition. 



A. The exact factorization 

Theorem I. (a) The exact solution ofEq. dlQl > can be writ- 
ten as a single product 

*(r 1, E g, t) = *g g(i 1' *)x(E g, t) (23) 

where $r „ (r s, t) satisfies the Partial Normalization Condi- 
tion (PNC), ~ 



E 



J dl|%z(£|^)| 2 = 1 , 



(24) 



for any fixed nuclear configuration, R q_, at any time t. 

The PNC is critical in making this theorem meaning- 
ful: Eq. j23l on its own would be rather meaningless, 
because, for example, one could then simply just take 
x(R <x, t) = 1. In fact, one can come up with many dif- 
ferent decompositions that satisfy Eq. d23l but that vio- 
late the PNC Eq. {24); it is the latter that makes the de- 
composition unique up to a gauge-like transformation, 
as we shall see shortly in Section III Bl We will also see 
there that it is the PNC that allows the interpretation of 
$R £ (r s, t) as a conditional probability amplitude, and 

x(Rff, t) as a marginal probability amplitude, leading to 
their identification as electronic and nuclear wavefunc- 
tions respectively. First, we prove Part(a) of Theorem I. 

Proof: Given ^(r s, R a, t), the exact solution of the full 
TDSE llT0t , we choose x(R <L, t) an d ®Ra(L s, t), at any 
instant in time, as 



X (R g, t) = e iS( l ^EJ#fe|lpP (25) 
and 

% £ (g|,t) = *(r|,S£,t)/x(Sg,t) (26) 

where S(R g, t) is real. The PNC Eq. ((24) then follows 
immediately: 



J «fe|*gg(£i>*)| ! 



EJdg\*(Lg,Rg,t) 

lx(R^)l 2 
lx(R<i,*)l 2 



lx(Eff,*)l s 



1. (27) 



This concludes the proof of Theorem I (a). It will be- 
come clear throughout this paper that, in many respects, 



4 



the nuclear factor x(E £, t) can be viewed as a proper 
nuclear wavefunction. Like in the static case [20], intro- 
ducing the phase factor in Eq. | |25)| allows x(E g_,t) to 
have the correct antisymmetry if the nuclear subsystem 



contains identical fermionic nuclei. 

Next comes the question; what equations do 
$Ro-(r s, t) and x(E £, t) satisfy? The answer entails 
the second part of Theorem I: 



Theorem I (b) The wavefunctions $R £ (rs, t) and x(E£j *) satisfy: 



(Hel(r§.,Qg,t) - <E £,*))% 5(1 = idt^gig^t), 



(28) 



a=l 



ifftere f/ze electronic Hamiltonian is 



Hd (r s, R £, t) = Hbo (r, R, t) + Kxt fe *) + 



coup 
en 



$Ra,X 



(29) 



(30) 



Here £//e electron-nuclear coupling potential U t 
are 



coup 
en 



4> 



R<r,X 



jjcoup 
^ en 



a — 1 

e(E g,t) = E (% I ^ (d & 1 £ *) _ 0i 

s 

A (B_£,i) = E (%sW -«V a $ gg (t) 



, scalar potential e(Rer, i), and vector potential A a (Ra, t) terms 



X(E£,*) 



$R £ (t) 



A Q (E£,i)J • (-*V a - A a (Eg,t))J(31) 

(32) 

(33) 



if/iere (..|..|..) r denotes an inner product over all spatial electronic variables only. 

I 

Proof In order to derive the equations of motion for a functional of the time-dependent wavefunction ^>(t) 

and its complex conjugate. The equation of motion of 



<1>rct(i s,t) and x(E <L>t) we follow the strategy em 
ployed in the static case (see ref. j20ll ), i.e. we plug the 



product ansatz in the variational principle and search requiring the variation of the action 5 with respect to all 
for the stationary point. Afterwards we prove: if 
$r 2 (r s, t) and x(R £) t) are the solutions of Eqs. (|28l 
and ll29l , then ^RctCt. s, i)x(E £>*) is tne solution of 



the quantum system, the TDSE of Eq. flOl l, is obtained by 
requiring the variation of the action S with respect to all 
wavefunctions that satisfy the boundary condition 

(35) 



TDSE (flOb . We begin the derivation by briefly reviewing 
Frenkel's stationary action principle as this is the key in- 
strument to derive the equations of motion for each sub- 
system. 

The quantum mechanical action is defined as 



5V(ti) = 5V(t f ) = , 
to be stationary, i.e., 



6y*S = 0. 



(36) 



(34) 



Now we apply this general variational principle to our 
problem in the following way. We insert the product 
wavefunction in the action functional |34j, with Hamil- 
tonian given by Eq. Q}, rewriting it as 



<s[% s ,$| g ,x,x*] = E / J d =J d = 



\x\ 2 $k g ( H BO + Kit + E ^mJ " ^ J *lg 



( E 2m" + ^ " idt J x + lxl 2$ |^E]frHV Q x/x) • HV„<l» a£ ) 



(37) 
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The equations of motion for $r „(r s,t) and x(R a, t) are obtained by requiring the action functional ll37t to be 
stationary with respect to variations of each wavefunction subject to the PNC l(24jl ; i.e., 

(38) 



Variation of Eq. 07t with respect to a (r s) leads to 
Ixl 2 [hbo + tit + E l5 " 04 ) 



^[*R a, $R a , X, X*] *5[*R 2, <&R CT , X, X*] 

= and ~- == ; = 



idt X 



<I> 



R a 



(E j^HVax/x) ■ H v «%g)l 



Dividing the expression above by |x| 2 and rearranging yields: 



H 



BO + Vext 



M, 



■ %g 

(39) 



Variation of Eq. | |3"7t with respect to x* yields 



E ^rHVax/x)- a « 



x = «ftx 



(40) 



where we enforced the PNC, and defined product $r „ (r s, t)x(R <r, t) is an exact solution of the 

I- TDSE. Approximate solutions of the TDSE may satisfy 

A q [$ro_] := E / ^i^R. o-(l|)( — *V q <I>r o-(lg)) (41) the stationary action principle, if variations are taken 
s J over a limited set of wavefunctions, e.g. the multi- 

configuration time-dependent Hartree equations HH 
may be derived via the Frenkel variational principle. 



This is a real- valued vector potential (see shortly). In- 



serting Eq. CP on the RHS of Eq. © leads, after some To dispel any possi ble doubts that the product form 
straightforward algebra, to Eqs. (28J33J. The product of Eq gjj) subject to Eq. ® is general, we now ver- 



wavefunction Eq. | |23l , satisfying these equations, there 
fore represents a stationary point of the action func- 
tional ( f37b under the PNC Eq. l l24l . To complete the 
proof, it remains to verify that if ^r^Cl £>^) satis- 
fies Eq. ll28l and x(R 2L,t) satisfies Eq. ll29l , then the 



ify that our solution is exact and not an approxima- 
tion. Applying the product rule, idt^(r s, R a,t) = 
x(R g, t)id t ^ g (L g, t) + $fig(E |> *) i9 *xTS iL , and 
inserting Eqs. 11281 and (|29|l , we obtain 



X N<9*% g J = x [Hbo + Kxt ) % g + X E 



HVq - A a ) S 

2M„ 



+xE 



(-»VaX/X + A a) • (~iVa ~ A a ) 



R , 



X£*R, 



(42) 



(-iVa + A a (Ka,t)) 2 
$ g ^ ( iaeX ) =$r^E 2AJ a X + % gT^ t X + % g^X 



(43) 



Summing Eqs. | |42|| and j43l leads to the TDSE for the complete system and completes the proof that the wave- 
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functions satisfying Eqs. (|28H33b do solve the TDSE ex- 
actly. 

Alternatively Eqs. J2.8H33|> can be obtained by replac- 
ing *(r || g,t), in the TDSE ((10) , by the product 
^Ra(|"|,*)F(Rg)*) and using the PNC d23J- The form 
of electron-nuclear coupling term, Eq. 1)3 It , is the same 
as the static case (see ref. l20ll ). The exact TDPES, 
Eq. d32l l, on the other hand is not simply the expecta- 
tion value of H e i but contains, in addition, the term 
< a_\ — i<9 t $R a_ >■ The appearance of this term is 
essential to ensure the form invariance of the Eqs. j28T 
[33b under the gauge transformation |(44jl that will be dis- 
cussed in Section IIB. 



B. Uniqueness of the electronic and nuclear 
wavefunctions 

We now delve a little deeper into features of our ex- 
act factorization. As briefly mentioned earlier, the fac- 
torization can be viewed in a standard probabilistic set- 
ting JH: The square of the molecular wavefunction 
can be viewed as a multivariate probability distribution, 
that can be factorized into a marginal probability of a 
set of variables (the nuclear coordinates) and a condi- 
tional probability of the rest of the variables (the elec- 
tronic coordinates, conditionally dependent on the nu- 
clear coordinates). In this sense we identify y(R £>t) 
as the nuclear wavefunction (marginal probability am- 
plitude), and $r. o- (r s, t) as the electronic wavefunction 
(conditional probability amplitude). An equivalent for- 
malism is to view, instead, the nuclear wavefunction as a 
conditional probability amplitude depending paramet- 
rically on the electronic coordinate, i.e. \r s (R <l, t), with 
the electronic wavefunction as the marginal probability 
amplitude of the electronic coordinates, i.e. $(r s, t). 
We choose to use the former decomposition however to 
later make natural connections with the BO approach. 
In this section we argue why we can view the proba- 
bility amplitudes x(R a, t) and ^R^d s, t) as nuclear 
and electronic wavefunctions, and we will assign some 
meaning to the terms that arise in their equations of mo- 
tion. 

A first question that arises is: is this decomposition 
unique? We answer this in Theorem 2. 

Theorem 2 (a) Eqs. (I28H33P are form-invariant up to 
within the gauge-like transformation: 

x(Rg,t) :=e- tf (S^ X (E£,i) (44) 

A a (Rg,t) -> A (Sg,t) = A Q|g,i) + V a 0(Jig;,t) 
e(R g, t) -> e(E g, t) = e(R g, t) + d t 0{R g, t) (45) 

(b) The wavefunctions ®Ra(L s, t) and x(R <L,t) are 
unique up to within the (R a, t)-dependent phase transfor- 
mation, Eq. dill . 



To prove part (a), simply substitute Eqs. j44|l and j45|l 
into Eqs l|2"8"ll-l|3"3"ll. Part (b) is readily shown by first as- 
suming that $r aX and $r „x are f wo different repre- 
sentations of the exact wave function $(rs,R a, t) i.e. 

tffc&Eg,*) = %2(y,*)x(i^*) = ®Rz(Lg,tMRg,t) 

(46) 

Then 

y $Ra 

- = ^r=:.9E£,( (47) 

X % 2 

and 

|%^(l|,t)| 2 - |5(ga,t)| 2 |*aff(rg,*)| a . (48) 

From Theorem 1, both $r ^(r s, i) and $r £ (r s, t) sat- 
isfy the PNC. Hence, ~~ ~ 

J ^l%^(ii,i)l 2 = |g(Egi)l 2 ^ / *l*s s (*i>*)l a 

s s 

(49) 

and |5(Rct, t)| 2 = 1. Therefore, p(R er, i) must be equal 
to a purely (R er, independence phase: 

g(Rg,t) =e iS( i?' t) . (50) 

This completes the proof of theorem 2. 

The interpretation of <&r and x as electronic and nu- 
clear wavefunctions follows from the following obser- 
vations. The probability density of finding the nuclear 
configuration R at time t, J2 S J\^(L s, R t)\ 2 dr = 

lx(R <Zj t)\ 2 > as can readily be shown by substituting the 
product wavefunction Eq. j23l into the left-hand-side 
and using the PNC Eq. C2}. Not only does x(R_ g,t) 
therefore yield the nuclear (iV„-body) probability den- 
sity, we shall see later in Section IIII Al that it also re- 
produces the exact nuclear (iV„-body) current-density. 
The modulus-square of the electronic wavefunction, 
l%^(l|,i)| 2 = \V(Lg,Rg,t)\ 2 /\x(Qg,t)\ 2 gives the 
conditional probability of finding the electrons at r with 
spin configuration s, given that the nuclear configura- 
tion is R a. 

Note that, strictly speaking, the definition of the 
conditional probability amplitude |3>Rff(r s, t)\ 2 via 
Eq. ll26l , only holds for non-zero marginal probabilities 
lx(R SLit)\ 2 - m the case the nuclear density, and the 
full molecular wavefunction, have a node at some R Q , 
the electronic wavefunction would be defined by tak- 
ing a limit. However, it is actually very unlikely that 
the nuclear density has a node 13d, l3lll . This can be 
seen by expanding the full electron-nuclear wavefunc- 
tion, ^(r s,R a, t), in terms of the BO-electronic states, 
as in Eq. \17\ . Then, the nuclear density can be expressed 
as an infinite sum of non-negative terms: 

oo 

lx(§£,i)| a = Elx?°(S£,*)| a - < 51) 

3=1 
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In general, it is extremely unlikely that every term in the 
summation becomes zero at the same nuclear configu- 
ration R Q a , unless dictated by symmetry |20|| (see end 
of this section for a discussion on symmetry). Symme- 
try dictated nodes likely lead to a finite, well-defined, 
value of |$r ff(r s, t)\ 2 due to the linear behavior of the 
wavef unctions in the vicinity of these nodes. 

Eqs. <f28t— d33l > determine the exact time-dependent 
molecular wavefunction, given an initial state. As writ- 
ten, the nuclear equation is particularly appealing as 
a Schrodinger equation with both scalar and vector- 
potential coupling terms contributing effective forces on 
the nuclei including any geometric phase effects. We 
call e(R g_, t) and A(R a, t) the exact TDPES and exact 
time-dependent Berry connection, respectively. These two 
quantities, along with the electron-nuclear coupling po- 
tential C/g° wp [$R a, x]/ mediate the coupling between the 
nuclear and the electronic degrees of freedom in a for- 
mally exact way. The three sections in Section |lll] are 
each devoted to a closer study of these terms. 

We conclude this section by discussing the symme- 
try properties of x(R a_, t) and "Jr a (r s, t): The nuclear 
wavefunction x(R a, t) must preserve the symmetry 
of the full electron-nuclear wavefunction >t(r s,R a, t) 
with respect to exchange of identical nuclei. This 
constrains the allowed gauge transformation (04)- 
105}. The electronic wavefunction &ro_{l t) = 
f(r s,R a, t) /x(R 3., t) is invariant under any nuclear 
permutation because any fermionic sign cancels out be- 
tween the full molecular wavefunction and the nuclear 
wavefunction. 

In the rest of the paper, we drop the spin indices g_ and 
s for notational simplicity. 



and n = M/(M + 1) is the reduced mass. The full wave- 
function, Eq. represents free-particle plane-wave 
motion in the center of mass coordinate, with K repre- 
senting the total momentum of the system. The form of 
Eq. (|53l l suggests one possible factorization for Eqs. j23|l 
-dtlj as: 

if -K 2 t i M ir Ti\ 
X(R,0 = e 4 UM + l) + (Af + l)"-"-J 

$R.(r,i) = e*wft)0(r - R,t) (55) 
with the exact Berry potential and TDPES given by 

A(R, t) = —i [ <f>* (r - R, t)V R (/)(r - R, t)dr = (56) 



e(R,*)= 2(M+1) + R-E(t). (57) 

The vector potential, Eq. f56l l, is zero in the gauge im- 
plicit in our choice for Eqs. ll55l l. This is easily confirmed 
by inserting Eqs. (55) in the nuclear equation d29l l, which 
reads for our problem, 

(-iV + A) 2 - R • E(i) + e(R, t)j X (R, t) = Wtxfr, *) 

(58) 

Eqs. | |57t and ( |5Hl l show that, in this case the role of the 
TDPES is to cancel out the external laser field in the nu- 
clear equation, which is exactly as it should be. Only 
by this cancellation the nuclear motion can be a plane 
wave. 



III. THE EXACT ELECTRON-NUCLEAR COUPLING 
TERMS 



C. Simple Illustration: the H atom in an electric field 

The example of the hydrogen atom in an electric field 
provides a simple demonstration of our formalism. The 
Hamiltonian is 



1 



1 2 1 2 

H = ~2M VR ~2 Vr ~ |R-r| 



(r-R)-E(t) (52) 



where r and R are the electron and proton coordinate 
respectively, E(t) is the applied electric field, and M is 
the proton mass. The exact solution is known: in terms 
of the center of mass and relative coordinates, Rcm = 
(r + MR) / (M + 1) , u = r — R, the problem is separable, 
and we have 



'J' (Rcm, u, t) = e H K ' RcM_ 2(17+1) 'J ^ t ) ( 53 ) 
where </>(u, t) satisfies the following equation: 



- - + u • E(t) 0(u, t) = id t 4>(n, t) (54) 



We now take a closer look at each of the three terms 
A(R g,t), e(R g,t), and E/^[$ B2 , x], that mediate 
the coupling between electron and nuclear dynamics ex- 
actly. In these three terms, all of the non-adiabatic cou- 
pling effects of the Born-Oppenheimer expansion are ef- 
fectively contained. 



A. The time-dependent Berry connection 

Eqs. ||2"Bll-(|55t demonstrate that a Berry connection in- 
deed appears in the exact treatment of coupled electron- 
ion dynamics, a question which was raised in the intro- 
duction. In this section, we point out some properties of 
this object to help us understand what it represents. 

First, we show that the vector potential A a is real. 
Taking the gradient with respect to nuclear coordinates 
of the PNC (Eq. ©), yields 







dr$ R (r)$ R (r) 



2Re J dr$|(|)V a $s(E) (59) 
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(using the product rule). Comparing with the definition 
Eq. (f33b , we readily conclude A Q is real. 

Second, we insert Eqs. (25l l and (26^ into Eqs. ((33b to 
reveal the following expression for the vector potential: 

/m(tf(t)| V a *(t)> P 
A «®> *) = |x(R f)|2 — E - V Q S(E, *) (60) 

This shows that the vector potential is the difference 
of paramagnetic nuclear velocity fields derived from 
the full and nuclear wavefunctions. In fact, since 
Im ($(t)\ V a ^{t)) r is the true nuclear (many-body) cur- 
rent density, Eq. d5Ut implies that the gauge-invariant 
current density, im(x*V a x) + |x| 2 A a , that follows from 
the nuclear Hamiltonian in Eq. (29) does indeed repro- 
duce the exact nuclear current density [32]. As dis- 
cussed in the previous section, the solution x(R>£) of 
Eq. ll2"8t yields a proper nuclear many-body wavefunc- 
tion: Its absolute-value squared gives the exact nuclear 
(A^-body) density while its phase yields the correct nu- 
clear (A rl -body) current density. (The nuclear kinetic en- 
ergy evaluated from x(R, t) does not equal the nuclear 
kinetic energy evaluated from the full molecular wave- 
function, and their difference is determined by U^° up , as 
will be discussed in Section PlIQ . 

Another interesting aspect of expression 1(601 is that it 
can help to shed light on the question of whether the 
exact Berry potential produces a real effect or whether 
it can actually be gauged away by a suitable choice of 
9(R,t) in Eqs. (45H55) . Provided the phase S(R,t) is 
spatially smooth, the last term on the right-hand-side of 
Eq. d60b can be gauged away so any true Berry connec- 
tion (that cannot be gauged away) must come from the 
first term. In the conventional analyses of conical inter- 
sections, the phase may not be smooth: for example, in 

I 



The exact Berry potential is thereby expressed as a lin- 
ear combination of the diagonal and off-diagonal BO 
derivative couplings. Any gauge-invariant part of the 
Berry connection, that would give rise to a non-zero 
Berry phase, arises from the part of Eq. [64] that has a 
non-zero curl. In the case of a real-valued electronic 
wavef unction, each of the three terms of Eq. [64] vanishes 
independently giving rise to a zero vector potential. 



the Herzberg and Longuet-Higgens model (26l, HH], the 
two (single-valued) nuclear wavefunctions associated 
with a two-state conical intersection between traditional 
BO surfaces, each have a phase S = ±0/2, undefined at 
the origin. This has a singular gradient, yielding a delta- 
function at the origin in the curl of the vector potential, 
thus contributing a non-zero Berry phase. Whether a 
similar effect occurs for the exact time-dependent nu- 
clear wavefunction remains to be explored. When the 
exact \&(i) is real- valued (e.g. for a non-current-carrying 
ground state) then the first term on the right-hand-side 
of Eq. 11601 vanishes and hence gives a vanishing contri- 
bution to the exact Berry connection. Whether, and un- 
der which conditions, the full Berry connection d60t can 
be gauged away remains an open question at this point. 

Finally, it is also instructive to express the vector po- 
tential in terms of the BO electronic basis states of Sec- 
tion llCI We first expand the electronic wavefunction: 

oo 

%(E,f)=^C 3 -(E,t)4(l) (61) 
i=i 

where orthonormality of the R (Eq. (16|| ) means 

Cj (R, t) = j d^g (r)*s(E. *) ■ (62) 
The PNC condition becomes 

oo 

J2\Cj(R,t)\ 2 = l (63) 
j=i 

Inserting Eq. [ |6"ll into Eq. J33b , and noting the definition 
of the non-adiabatic derivative couplings of Eq. [22] we 
obtain 



(64) 



I 

B. The Time-Dependent Potential Energy Surface 

The time-dependent potential energy surface e(R, t) 
of Eq. (32|l provides an exact time-dependent general- 
ization of the adiabatic BO potential energy surface. As 
such, it should prove to be a powerful interpretive tool 
for general time-dependent problems. This will be ex- 
plored in section IV. We now begin by analyzing the ex- 
pression Eq. ( T52l in a little more detail. 

First, consider the expectation value of the electron- 
nuclear coupling term, ($R|C/ e c °" p |$ a ) of Eq. (3D that 
appears in the TDPES. Only the first term of Eq. (3D con- 
tributes to the expectation value: the second term goes 
to zero, due to the very last parenthesis, ($r| — iV a — 
A a (R, f)|$R), which vanishes due to the definition of 
the vector potential. So we have 



A Q (R, t) = J2 ( - *Cj GL> *) V«C, (R, t) + | Cj (R, t) \ 2 ^° a (R) + ]T C* (R, t)C (R, t)F%% (R) 



m,t) = [{^\H BO +v: xt (r,t)\^ m ) r= -i{^\d t ^) g +j^ 

= (<%|#bo + v: xt (*, t)|*a) £ - »<%lfc%> £ + E 



($RlHV Q -A Q (E,i)) 2 |3>R> 



(v a %|v„*s) 

2M~ a 



E 



2M Q 



(65) 



where the second line results from expanding the square an expression for the TDPES as an expansion over BO 

in the first, and making use of the definition of the vector states. Inserting Eq. i EjTT i into Eq. [65] and performing a 

potential. little straightforward algebra, we obtain 
As we did for the vector potential, we now provide 



e(E, i) = E |Q(E, t)\ 2 Vi (R) + £ C*(E, i)C*(E, t)(^|K e ^fc *)|4) S - E iC 7d< *) 

+ E E I V° C i I" + E ( iV « • T ™ ~ - 2 E Irn( C iV*C*F?° a ) Al (E, «) I (66) 

a ™ \ i ji ji J 



(the expansion of the last term A 2 may be obtained from 
Eq. [64)l . Notice that all the BO surfaces, as well as non- 
adiabatic couplings, are contained in the the exact TD- 
PES. 



where 



(69) 



C. Electron-Nuclear Correlation 

The TDPES and Berry connection discussed in the 
previous two sections directly determine the evolution 
of the nuclear wavefunction (Eq. l|29l l), containing the 
effect of coupling to the electrons in an exact way. 
The electron-nuclear coupling term E/£n™ p enters the nu- 
clear equation indirectly via its role in determining <f>R 

through Eq. (gg) and (SO). Eq. C2) expresses U™ up as a 
functional of the electronic and nuclear wavefunctions, 
and now we shall derive another expression for it that 
shows that it measures the difference between the nu- 
clear kinetic energy evaluated from the full wavefunc- 
tion and that evaluated on the nuclear wavefunction. 
We isolate the term involving U^° up in Eq. (28) , and in- 
sert $r =^/x- This leads to: 



R 



R 



■K e xt +e(R,t) (67) 



Next we insert in Eq. f&7$ the TDSE ^ and Eq.®, sat- 
isfied by ^ and \ to obtain 



jjcoup 



$R,X 



$R(r,i) f y f 

— — ± n^ ± n 



(68) 



Multiplying Eq. d68l by I^r^IxI 2 and integrating over 
all coordinates leads to: 

(*|r„|*)r 1 R-(x|r n |x>a = / ^lx(E,*)l 2 <$al^T p l<&a} 



(70) 

This means the nuclear kinetic energy evaluated from 
the full molecular wavefunction, and that evaluated via 
the expectation value of the nuclear kinetic energy op- 
erator in Eq. 1 1291 on the nuclear wavefunction are not 
equal: their difference is given by the nuclear-density- 
weighted integral of the electron-nuclear coupling po- 
tential. 



IV. MODEL OF H+ IN A LASER FIELD 

In this section, we illustrate the usefulness of the TD- 
PES using a simple, numerically exactly solvable model: 
the H% molecular ion subject to a linearly polarized 
laser field. By restricting the motion of the nuclei and 
the electron to the direction of the polarization axis of 
the laser field , the problem can be modeled with a ID 
Hamiltonian featuring "soft-Coulomb" interactions lf33l — 



H(t) 



i d 2 



i d 2 

MdR 2 ~ ^dx 2 
1 



1 



V0.03 + R 2 
1 



+ V l {x,t) 
(71) 



y/l + (x- R/2) 2 y/1 + (x + R/2y 
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where R and x are the internuclear distance and the elec- 
tronic coordinate as measured from the nuclear center- 
of-mass, respectively, and the electronic reduced mass 
is given by \i e = (2M)/(2M + 1), M being the pro- 
ton mass. The laser field is represented by Vi(x,t) = 
q e xE(t) where E(t) denotes the electric field amplitude 
and the reduced charge q e = (2M + 2)/(2M + 1). One- 
dimensional soft-Coulomb atoms and molecules have 
proven extremely useful in the study of strong-field 
dynamics since they allow numerically accurate solu- 
tions to problems involving correlated electron dynam- 
ics as well as correlated electron-nuclear dynamics that 
would be computationally far more demanding for the 
full three-dimensional atoms and molecules, while cap- 
turing the essential physics of the latter, e.g. multi- 
photon ionization, above-threshold ionization and dis- 
sociation, enhanced ionization, non-sequential double- 
ionization, high-harmonic generation, and non-BO ef- 
fects (e.g. Refs. IZl l35l - l41|| '). We study the dynamics of 
the model # 2 + system under a A = 228 nm (5.4eV) UV- 
laser pulse which is represented by 



E(t)=E f(t)an{ut), 



(72) 



with two peak intensities, h = \Eq\ 2 = 10 14 M / /cm 2 and 
h = \E \ 2 = 2.5 x 10 13 VK/cm 2 . With this frequency an 
energy that is about twice as much as the dissociation 
energy of the model molecule (2.8782eV) is achieved, so 
dissociation is expected. The envelope function f(t) is 
chosen such that the field is linearly ramped from zero to 
its maximum strength at t = T ramp and thereafter held 
constant (Fig. [T]l: 



/(*) 



t f^ramp 

l 



< t < T ra 

Tramp ^ t <C Ti Q f 



nip 



(73) 



The rise-time was chosen as T ramp = lOr while the total 
simulation time was T to t = 25t, where r = — denotes 



the optical cycle. 



0.06 




10 15 
Optical Cycle (0.76 fs) 



25 



FIG. 1. A = 228 run laser field, represented by E(t) = 
Eof(t) sin(ut), for two peak intensities, I\ = \Eo\ = 
10 14 VK/cm 2 and I 2 = \E a f = 2.5 x W 13 W/cm 2 . The envelope 
function f(t) is chosen such that the field is linearly ramped 
from zero to its maximum strength at t — 7.6 fs and thereafter 
held constant. The highlighted area represents the optical cy- 
cle that will be focussed on in later graphs. 

The same system and parameters were studied in 
Ref. 13^1 where the importance of electron-nuclear cor- 
relation was highlighted: a two-configuration corre- 
lated ansatz for the time-dependent electron-nuclear 



wavefunction was able to describe photodissociation 
processes in many cases, while a simple uncorrelated 
Hartree product of an electronic and a nuclear wave- 
function almost always failed. In the present work we 
analyse the dynamics via the numerically exact TDPES, 
finding it very useful in understanding and interpreting 
the motion. We note that the laser-field does not couple 
directly to the nuclear relative coordinate R, but only in- 
directly via the TDPES. 

Starting from the exact ground-state as initial con- 
dition, we propagate the TDSE numerically, using the 
second-order split-operator method [42], to obtain the 
full molecular wavefunction fy(x,R,t). As there is 
only one nuclear degree of freedom (after separating 
off the center-of-mass motion), we can fix the gauge in 
Eqs. <f44b— <f45l> such that the vector potential l l60T l vanishes 
identically. For one-dimensional problems this is always 
possible with the choice: 



Imfdx^*(x,R,t) d * { *£' t) 



(74) 



So we can calculate S(R,t), the phase of the nuclear 
wavefunction, as well as \x(R,t)\ 2 , the nuclear den- 
sity, from the computed exact time-dependent molec- 
ular wavefunction. Being equipped with the nuclear 
wave-function, x(R^) (= lx(-R> t)\e iSI - R '^) , we then 
compute the TDPES by inverting the nuclear equation 
of motion j29l l. 

We will compare the exact dynamics with the fol- 
lowing three approximations: (i) the usual Ehrenfest 
approximation, where the nuclei are treated via clas- 
sical dynamics, evolving under the force — VVe/ ir = 
-VsW n „(R) - / dgn(r, t)V^W en (r, ||), with n(r, t) be- 
ing the one-body electron density , (ii) the "exact- 
Ehrenfest" approximation, which substitutes the exact 
TDPES for the Ehrenfest potential VEhr in the usual 
Ehrenfest approach and, (iii) an uncorrelated approach, 
the time-dependent Hartree (self-consistent field) ap- 
proximation, ^fffoR, t) = 0(r, t)x(R, i), where the 
electronic part does not depend on R at all. This 
includes a quantum treatment of the nuclei, but no 
electron-nuclear correlation. 



1. High intensity: h = l0 14 W/cm 2 

The exact TDPES, along with the corresponding nu- 
clear density, \x(R,t)\ 2 , are plotted in Fig.|2]at six snap- 
shots of time. The initial TDPES lies practically on top of 
the ground-state BO surface, plotted in all the snapshots 
for comparison. 

The dissociation of the molecule is dramatically re- 
flected in the exact TDPES, whose well flattens out, 
causing the nuclear density to spill to larger separa- 
tions. Importantly, the tail of the TDPES alternately falls 
sharply and returns in correspondence with the field, 
letting the density out; the TDPES is the only potential 
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acting on the nuclear system and transfers energy from 
the accelerated electron to the nuclei. 
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FIG. 2. Snapshots of the TDPES (blue solid lines) and nu- 
clear density (black solid lines) at times indicated, for the 
H^" molecule subject to the laser-field with the peak intensity 
h — 10 14 W/cm 2 . The solid circles indicate the position and 
energy of the classical particle in the exact-Ehrenfest calcula- 
tion. For reference, the ground-state BO surface (red dashed 
lines) is shown. 



In Figure|3]we focus on six equally-spaced time snap- 
shots during the optical cycle shaded in Figure [T] The 
lower panel shows the TDPES, with its characteristic os- 
cillations, along with the nuclear density as a function of 
the internuclear coordinate, \x(R,t)\ 2 . The upper panel 
shows a color map of the conditional electronic proba- 
bility density, t)\ 2 , i.e. the probability of finding 
an electron at x at a fixed nuclear separation R. While 
at small internuclear distances (around and below the 
equilibrium separation) the electron remains localized 
in the middle between the two nuclei, at larger separa- 
tions one clearly sees the preferential localization of the 
electron density near the two nuclei, i.e. on one side or 
the other. At even larger separations we see streaks of 
ionizing electron density in both directions. For the full 
story, we must multiply the conditional probability den- 
sity of the upper panels with the nuclear density shown 
in the lower panel, to obtain the total electron-nuclear 
density; this is shown in Figure |4j indicating the prob- 
ability of finding, at the time indicated, an electron at 
position x and the nuclear separation R. 

The top left-hand panel of Fig. [5] shows the expecta- 
tion value of the internuclear distance 



<i?>= (*(t)|ii|tt(t)) 



(75) 



along with the results from the three approximate meth- 
ods described earlier. The lower left-hand panel shows 
the ionization probabilities. In principle, the latter re- 
quires projections of the full wavefunction on all contin- 
uum states which, in practice, are difficult to calculate. 



FIG. 3. Snapshots of the TDPES (blue lines), nuclear density 
(black) and the electronic conditional-density (color map) at 
times indicated during an optical cycle, for the molecule 
subject to the laser-field with the peak intensity h = 
10 14 W/cm 2 . For reference, the ground-state BO surface is 
shown as the red line. 
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FIG. 4. Snapshots of the total electron-nuclear density at times 
indicated during an optical cycle, for the molecule subject 
to the laser-field with the peak intensity h = 10 14 W/cm 2 . 



Alternatively, we use a geometrical concept 14711 , accord- 
ing to which the total ionization probabilities can be ob- 
tained from 



Pion(t) = 1 



box e 



dx / dR\^{t)\ : 



(76) 



The electrons leaving the "electronic analyzing box" 
(box e ) are thereby identified with ionized electrons. The 
ionization box here was chosen to be |x| < 10. The inter- 
nuclear distance together with the ionization probabil- 
ity support a Coulomb-explosion interpretation of the 
dissociation: first, the system begins to ionize, then the 
nuclei begin to rapidly move apart under their mutual 
Coulomb repulsion increasingly sensed due to weaker 
screening by the reduced electron density. Turning now 
to the approximations, we observe that all the methods 
yield dissociation and some ionization. The expectation 
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value of the internuclear distance in Fig.[5j demonstrates 
that among all the approximate calculations employed 
here, the exact-Ehrenfest is most accurate. Referring 
back to Figure |2 the solid circles indicate the classical 
nuclear position and energy of a particle driven by the 
exact-Ehrenfest force. One can see that it rapidly picks 
up kinetic energy above the TDPES, supporting the fact 
that the nuclear dissociation mechanism is an essentially 
classical one in this case. The exact-Ehrenfest calculation 
even does better than TD-Hartree which treats the pro- 
tons quantum mechanically thus showing the overar- 
ching importance of electron-nuclear correlation in this 
case. 
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FIG. 5. Dissociation and ionization for intensity I\ (left) and I2 
(right). Top panels: the internuclear separation (R)(t). Lower 
panels: The ionization probability. 

In fact, the Hartree description is worse than it may 
seem from just looking at the internuclear separation in 
Fig- El In Figure [6] we plot the time-dependent Hartree 
potential energy surface and Hartree nuclear-density. 
Both are dramatically different from the exact TDPES 
and exact nuclear density of Figure|2] At the initial time, 
the Hartree potential is reasonably good near equilib- 
rium but poor at large separations 13711 : this is a conse- 
quence of the conditional electron probability being in- 
dependent of the nuclear coordinate, and therefore only 
yielding a realistic result where the energy is optimized, 
which is at equilibrium separation. As time evolves the 
minimum of the Hartree surface moves out and begins 
to widen, cradling the nuclear density, which more or 
less retains its Gaussian shape, unlike the exact density; 
only at larger times does the surface open out. 

2. Lower intensity: I 2 = 2.5 x 10 13 W/cm 2 

We now consider the dynamics under a field of 
weaker intensity. Figure [7] plots the TDPES, whose tail 
displays similar oscillations as in the higher intensity 
case. The nuclear density appears to leak out to larger 
separations, although more slowly than in the previ- 
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FIG. 6. Snapshots of the time-dependent Hartree nuclear- 
potential (blue lines) and nuclear density (black) at times in- 
dicated, for the molecule subject to the laser-field with the 
peak intensity h — 10 14 W/cm 2 . For reference, the ground- 
state BO surface is shown as the red line. 



ous case; indeed from the right panels in Fig. |5j we see 
that the exact calculation leads to dissociation. How- 
ever, Fig. [5] (upper right panel)also shows that none of 
the approximations dissociate, in contrast to the previ- 
ous case. The Hartree and Ehrenfest methods also show 
negligible ionization, compared to the exact case; but 
even in the exact case the ionization probability is very 
small, indicating a different mechanism of dissociation 
than in the stronger field case. It may be at first surpris- 
ing that the exact-Ehrenfest calculation does not dissoci- 
ate the molecule, given that it is based on the exact TD- 
PES, however an examination of classical dynamics in 
the TDPES of Fig. [2] can explain what is happening. The 
solid dot in Fig.|2]indicates the classical position and en- 
ergy, and we see that it is always trapped inside a well 
in the TDPES, that remains at all times. This suggests 
that tunneling is the leading mechanism for the disso- 
ciation: a classical particle can only oscillate inside the 
well, while a quantum particle may tunnel out, as in- 
deed reflected in Fig. [5] Although the tail has similar 
oscillations as for 1\, this does not lead to dissociation 
of classical nuclei due to the barrier; the TDPES in this 
case transfers the field energy to the nuclei via tunnel- 
ing. Although the exact-Ehrenfest calculation shows a 
larger amplitude of oscillation than the others, it ulti- 
mately cannot tunnel through the barrier. 

As in the previous case, we plot in the top panels of 
Fig. [5] the electronic conditional density \Qr(x, t)\ 2 over 
one optical cycle, while the lower panels illustrate again 
the opening and closing of the TDPES as the field oscil- 
lates. Like in the previous case, for small R near equi- 
librium, the electron density is localized in between the 
nuclei, while for larger R, there is some polarization 
towards one side or the other. To get the full picture, 
one must multiply the top panels by the nuclear density 
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\x(R,t)\ 2 , to obtain the total electron-nuclear probabil- 
ity density, shown in Figure [9] It is evident in this graph 
that there is much less ionization than in the previous 
case, and the dissociation is slower. 
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Although the Hartree approximation treats the nuclei 
quantum mechanically, and therefore allowing tunnel- 
ing in principle, tunneling and dissociation do not ac- 
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FIG. 9. Snapshots of the total electron-nuclear density at times 
indicated during an optical cycle, for the molecule subject 
to the laser-field with the peak h = 2.5 x 10 13 W/cm 2 . 



tually occur. The reason for this is clear from the shape 
of the Hartree potential, plotted in Fig. [TOj the Hartree 
potential essentially retains its initial shape at all times, 
making very small oscillations near the equilibrium sep- 
aration. As in the more intense field case, this is due 
to its uncorrelated treatment of the electron-nuclear sys- 
tem: the electronic wavefunction at any nuclear config- 
uration is always the same, and is best at equilibrium 
since initially it is determined by energy-optimization, 
from where it does not deviate far, due to the weak field 
strength. Unlike in the stronger field case, the Hartree 
surface never opens out. Dissociation via tunneling re- 
quires both a quantum mechanical description of the nu- 
clei and an adequate accounting of electron-nuclear cor- 
relation. 

We do not expect the TDPES to be so different from 
the BO surfaces in all cases. For example, in the case 
of field-free vibrational dynamics of the H2+ molecule, 
where we start with a nuclear wavepacket displaced 
from equilibrium on the ground BO surface, we find 
the TDPES follows closely the BO surface throughout. 
The non-adiabatic couplings are weak in this case. The 
TDPES for field-free dynamics in other systems with 
stronger non-adiabatic couplings will be published else- 
where nu. 

The purpose of comparing the exact results with these 
methods (TD-Hartree, Ehrenfest and exact-Ehrenfest) 
was primarily to support the conclusions drawn from 
the exact TDPES regarding the dissociation mecha- 
nisms. An interesting question is how well do the 
more accurate approximate PES's proposed recently 
(e.g. Ref Irl6ll ) compare with the exact TDPES; this will 
be investigated in the future. 



V. CONCLUSIONS 

In this paper, we have shown that there exists a rig- 
orous factorization of the exact molecular wavefunction 
into a nuclear wavefunction and electronic wavefunc- 
tion, each of which retains the usual probabilistic mean- 
ing. The exact nuclear iV^-body density is |x(IL0| 2 
while I^R^r, t)\ 2 represents the conditional probability 
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FIG. 10. Snapshots of the time-dependent Hartree nuclear- 
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dicated, for the molecule subject to the laser-field with the 
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reference, the ground-state BO surface is shown as the dashed 
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of finding the electrons at r, given the nuclear configu- 
ration R. Equations ll2"8l — l|33l l are the equations of mo- 
tion that the electronic wavefunction and nuclear wave- 
function satisfy, and show explicitly how the electronic 
and nuclear systems are exactly coupled. These equa- 
tions enable the time-dependent potential energy sur- 
face (Eq. (|32l l) and the time-dependent Berry connection 
(Eq. | [53l ) to be defined as rigorous concepts, and we 
have discussed some general properties of them, and of 
the electron-nuclear coupling operator Eq. yl) . 

The example of the one-dimensional Hj" molecule 
in an oscillating electric field, solved numerically ac- 
curately, demonstrated that the TDPES is a powerful 
tool to analyze and interpret different types of dissoci- 
ation processes. By studying the shape and evolution 
of the TDPES, comparing classical dynamics in this ex- 
act potential to the exact quantum dynamics, we were 
able to distinguish whether the dissociation proceeded 
via nuclear tunneling or more directly in Coulomb- 
explosion. For this example, the TDPES is the only po- 
tential determining the nuclear dynamics, exactly con- 
taining the coupling with electronic dynamics. The ex- 
ample demonstrated the importance of capturing both 
quantum effects in nuclear motion and electron-nuclear 
coupling; the Hartree approach, for example, despite 
treating the nuclei quantum mechanically, was unable to 



capture dissociation via tunneling as the shape of its po- 
tential surface was completely wrong. Thus, for exactly 
solvable systems, the TDPES, and in more general cases 
than the one studied here, the geometric phase, can be 
very useful interpretative tools for dynamics. The calcu- 
lation of a TDPES has quite some history in the strong- 
field community, and several possible definitions of TD- 
PES have been proposed in the literature. The crucial 
point of our work is that it provides a unique defini- 
tion of TDPES (unique up to within a gauge transfor- 
mation): If one wants the TD many-body Schrodinger 
equation l l29l l to give the correct N-body density and cur- 
rent density of the nuclei, then the scalar potential and 
the vector potential must be given by eq. 1 1321 and <f33b - 
There is no choice apart from the gauge. That means that 
with any advanced technique that yields the TD molecu- 
lar wavefunction ^(r, R, t) one can evaluate the TDPES 
and Berry potential by first calculating the factors from 
Eqs. l|25]l-d26t and then evaluating the TDPES and Berry 
potential from Eqs. 



From a practical point of view, Eqs. d28t-||331 are not 
easier to solve than the time-dependent Schrodinger 
equation for the full electron-nuclear system. Rather 
they form the rigorous starting point for making ap- 
proximations, especially for the systematic develop- 
ment of semiclassical approximations. In the large- 
nuclear mass limit, the electronic equation reduces to 
Cederbaum's time-dependent BO approximation JlLll8||. 
Taking the classical limit for the nuclei in the large-mass 
limit, one retrieves the Ehrenfest equations with Berry 
potential JH] (see also ji^|4§]). Treating the nuclei clas- 
sically but retaining their finite mass, one finds correc- 
tions to the Ehrenfest equations that better account for 
non-adiabatic transitions 15011 . A direction for future re- 
search is to capture some nuclear quantum effects by a 
semiclassical or quasiclassical procedure |5ll HU], built 
on the exact foundational equations presented here. An- 
other direction would be to use the formalism as a pos- 
sible starting point to develop electron-nuclear correla- 
tion functionals in a density-functionalized version of 
the electron-nuclear problem l46l . A promising route is 
to develop a time-dependent generalisation of the opti- 
mized effective potential scheme proposed in | 
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